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Introduction : 

Most  of  our  work  in  19  81  depends  on  the  observation 
that,  if  processing  time  is  not  important,  better  statistical 
estimates  can  be  made  by  minimizing  a  least  squares  type  objec- 
tive function  by  trying  all  possibilities,  rather  than  by  employ- 
ing non-linear  transformations  to  permit  the  minimization  to  be 
done  analytically.   The  first  two  sections  document  the  extent 
of  the  improvement  in  various  circumstances.   The  improvement  is 
equivalent  to  several  decibels  gain  in  signal-to-noise  ratio. 

Section  three  documents  an  idea  for  determining  direc- 
tion from  three  measurements.   The  procedure  does  not  require 
initial  estimates  of  speed,  range,  or  signal  strength. 

Section  four  is  unrelated  to  the  first  three  sections. 
It  uses  an  equivalent  sweep  width  idea  to  emphasize  the 
difference  between  area  search  and  barrier  patrol. 


1.    Performance  of  various  least  squares  objectives 

1 . 1   Background 

Consider  a  situation  where  a  sequence  of  measurements 
C, /  Cy , . . .  of  some  phenomenon  are  made,  each  of  which  depends 
in  theory  on  a  set  of  parameters   p   as  well  as  some  known 
independent  variables.   If  the  parameters   p   are  unknown,  it 
is  natural  to  estimate  them  by  adopting  the  ones  that  make 
theoretical  and  actual  measurements  agree  as  closely  as  pos- 
sible.  The  resulting  estimates,  as  well  as  the  difficulty  of 
obtaining  them,  depend  on  what  is  meant  by  "closely".   The 
purpose  of  this  section  is  to  explore  the  tradeoff  between 
accuracy  and  difficulty  that  results  from  different  definitions 
of  the  word  in  a  particular  estimation  problem.   Specifically, 
we  assume  that  the  theoretical   ith   measurement  is 


(1)  C.  =  K  x  +  vt.  -  x. I    ,  where 


v  and   x   are  the  speed  and  initial  distance  to  a  target, 

x.   and  t.   are  the  position  and  time  at  which  the  sensor 

obtains  the   ith  measurement,  and  where   K   and   a 

are  two  positive  parameters  that  will  be  assumed  known  initially 

So   p  =  (x,v)   and   x.   and   t.   are  the  independent  variables. 

1.2   Direct  Search  Method 

A  natural  way  to  estimate   p  would  be  to  minimize  the 
Euclidean  distance  between  (C, , . . . ,C  )  and  (C, , . . . ,C  ) ,  where 


is  the  number  of  measurements.   So  we  define 


n        -   2 
(2)         E  (x,v)  =  I       (C.  -  C.)   ,  and  let 

n        i=l    1     x 


(x  TR,  v nTR)  be  the  point  that  minimizes  E  (x,v) .   If  the  only 
difference  between   C.   and   C.   is  one  of  additive  white 
noise,  then   x     and   v  T    have  the  virtues  of  being  maxi- 
mum  likelihood  estimates  of   x   and  v.   Since  the  function 
E  (x,v)  has  no  easily  exploitable  analytic  properties  other 
than  being  a  sum,   x      and   v      are  most  easily  obtained 
by  direct  search  of  some  set   S   (typically  a  grid)  of  values 
for  (x,v) .   With  each  measurement,  the  stored  values   E  (x,v)   ar< 
updated  by  employing  the  equation 


(3)       En+1(x,v)  =  En(x,v)  +  (Cn+1  -  Cn+1)2 


and  then   xnTD   and   v_T_.   are  selected  by  finding  the  smallest 
stored  value.   The  disadvantage  of  this  method  is  evidently 
one  of  time,  particularly  if   S   is  a  large  set  (fine  grid) . 
The  advantages  are  that   E  (x,v)   is  essentially  the  right  thing 

to  minimize  when  additive  noise  is  the  dominant  difference 
between   C.   and   C.  ,  and  that  the  amount  of  time  required 
per  estimate  is  independent  of  the  number  of  measurements  that 
have  been  made. 


1.3   Least  Squares  Estimates 


Let     y.    =  C-1/"    '      h  •  ~cll/a  m  K"Valx  +  vti   "   xii    '    and 

J  2.  1 


n 
(4)  F  (x,v)  -J   (y   -  y  )2  . 

n        i=l    1     1 


In  circumstances  where  the  sign  of   x  +  vt-x.   is  known  and 
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constant  for  all  observations,   y.   is  a  linear  function  of 
x   and   v,  and  F  (x,v)  can  be  minimized  analytically.   In 
fact,  the  minimizing  point  (x  _  ,  v    )  can  be  obtained  using 
standard  linear  regression  formulas.   We  will  refer  to   x, qn 
and   v  _    as  "least  squares  estimates",  even  though   x^... 
and   vnTD   are  equally  deserving  of  the  name.   The  advantage 

1J1  IS. 

of  the  least  squares  estimates  is  the  simplicity  of  obtaining 
them.   The  disadvantages  are  that  they  are  not  maximum  likeli- 
hood if  it  is  the  quantity   C.   (rather  than  y.)  that  suffers 
from  additive  noise,  and  that  the  target's  direction  of  travel 
must  be  known. 


1.4   A  Comparison  of  Least  Squares  and  Direct  Search  Estimates 

A  simulation  was  performed  to  compare  (x    ,  v    ) 
with  (xnTD/  v^TTJ  .   The  true  (x,v)  was  (30NM,  5KT)  in  all 

1J1 K     Ul  i\ 

cases,  with   K  =  a   =    1.   The  independent  variables  were 

t.  =  i/6   and   (xi ,  x?,...)  =  (2,2,4,4,6,...),  which  corres- 

ponds  to  jerky  sensor  motion  at  an  average  speed  of  6KT.   The 

measured   C .   were  obtained  from  C .   by 
i  i    J 


(5)  Ci  =  |Ci  expCA^)  +  A2W2|  , 

where  W,   and  W-   are  unit  normal  samples  and  A,   and  A~ 
are  parameters.   The  additive  noise  is  governed  by   A?  ,  and 
the  multiplicative  noise  by  A,  .   The  reason  for  introducing 
multiplicative  noise  was  to  test  robustness  of  the  two  esti- 
mators to  the  form  of  the  noise. 

The  DIR  procedure  was  programmed  to  examine  only  the 
possibilities  1,2,...,  50  for   x   and   . 5 , 1 . 0 , . . . , 15  for   v  - 
a  total  of  1500  possibilities  -  from  which  the  best  was 
selected.   There  is  no  reason  in  principle  why  the  DIR  pro- 
cedure should  not  also  search  negative  values  and  thereby 
determine  the  sign  as  well  as  the  magnitude  of   x   and   v  . 
In  fact,  the  original  version  of  the  procedure  did  exactly 
that.   Initial  runs  with  this  version  revealed  that  it  some- 
times got  the  sign  wrong,  particularly  with  only  2  or  3 
measurements.   While  this  was  not  surprising,  it  was  surpris- 
ing (initially)  that   X     did  not  suffer  similarly.   But 
the  LSQ  procedure  does  not  mistake  the  sign  because  an  assump- 
tion about  sign  is  intrinsic  to  the  procedure,  and  the  assump- 
tion happens  to  be  correct  in  the  simulation;  the  LSQ  proce- 
dure showed  no  capability  for  estimating  (x,v)  when  that 
assumption  was  wrong.   The  "fair"  thing  to  do  seemed  to  be  to 
give  the  correct  direction  to  both  procedures,  which  is  why 
DIR  searched  only  positive  values  for   x   and   v  even  though 
this  would  not  be  true  in  an  operational  version. 


It  was  also  found  that  the  LSQ  procedure  occasionally 
produced  wild  estimates  that  would  simply  be  dismissed  in 
practice,  but  which  significantly  affected  the  summary  sta- 
tistics.  The  remedy  for  this  was  to  make   0   estimates  if 

XLSO  or  VLSO  exceeded  10°  or  50  i-n  absolute  value. 

Figures  1-16  show  the  results  of  comparing  the  LSQ 
and  DIR  procedures.   Each  figure  displays  the  experimental 
mean  and  standard  deviation  (in  the  +/-  column)  of  the  four 
estimates   X  .  ,  V    ,  X    ,  and  V  T  ,  on  each  of  20  successive 
measurements.   All  statistics  are  based  on  50  independent 
replications,  so  a  68%  confidence  interval  would  have  a  width 
of  whatever  is  in  the  +/-  column  divided  by  /50  ~  7. 

Figures  1-8  show  performance  with  no  multiplicative 
error  and  increasing  amounts  of  additive  error.   The  scale  of 
A2  was  chosen  to  cover  the  range  from  good  performance  to 
disaster;  with  different  values  for   v   and   x,  a  different 
scale  would  be  needed.   When   A2  =  .002,  both  procedures  obtain 
good  estimates  on  the  average.   Note,  however,  that  the  LSQ 
procedure  has  a  higher  standard  deviation.   On  the  second 
measurement,  for  example,  X  c    has  a  standard  deviation  of 
3.98,  which  considerably  exceeds  1.88  for   X_T_..   Neither  pro- 
cedure  is  of  much  use  for  estimating  velocity  after  only  two 
measurements.   As  A2  is  increased,  the  fact  that  LSQ  tends  to 
have  larger  standard  deviations  than  DIR  is  accentuated.   By 
the  time   A2  =  .02,  LSQ  is  probably  useless,  whereas  DIR  might 
still  have  some  value  for  estimating   x   or  even   v   with 
enough  observations. 


Figures  9-12  show  results  with  A2  =  0  and  increasing 
values  of  Al .   Again,  LSQ  seems  to  have  greater  standard 
deviations  than  DIR. 

Figures  13  and  14  show  results  with  mixed  error  types. 

Figures  15  and  16  show  results  where  the  true   X   is 
30.5,  which  doesn't  happen  to  coincide  with  any  point  that 
DIR  searches.   In  figure  15,   X      has  a  larger  standard 
deviation  that   xLq_,  presumably  because  DIR  is  always  at 
least  .5  away  from  the  truth.   Figure  16  shows  that  this  effect 
is  negligible  when  measurement  errors  are  larger. 
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1.5  The  case  where   K   is  unknown 

All  results  reported  so  far  have  been  for  the  case 
where   K   and   a   are  known.   In  Sections  1.5  -  1.7,  we  ex- 
plore the  consequences  of  assuming  that   K,  x,   and   v   must 
all  be  estimated,  with  only   a   being  known.   We  should  anti- 
cipate that  good  estimates  of   x   and   v   will  now  require 
less  noisy  measurements  than  before,  since  there  is  more  esti- 
mating to  be  done.   The  principal  object  remains  the  comparison 
of  the  direct  search  and  least  squares  procedures. 

1.6  Direct  Search  Method  (K   unknown) 

Rewrite  (1)  as   C  =  Kd .    to  emphasize  that   C-   is 

11         ^  l 

a  linear  function  of   K  .   So   d.  =  |x  +  vt .  -  x,     .   The 
object  now  is  still  to  minimize   E  (x,v)   as  given  by  (2) ,  but 
in  the  current  notation 


n  ?  2 

(6)     E  (x,v)  =       I        (Ci  -  Kdi)   =  An  -  2KBn  +  K  Dn  , 

i=l 


n   ,2 


where        A   =  j      C .  , 


n 
B   =   7   C.d.  ,  and 
n    i=l   X  X 


n    2 


d  =  I    a* 
n      iii    x 
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Even  though  it  would  be  possible  to  minimize   E  (x,v)   by 
setting  up  a  three  dimensional  grid  on   K,  x,   and   v  ,  it  is 
easier  to  do  the  optimization  on  K   analytically,  since  (6) 
is  simply  a  quadratic  expression  in   K  .   By  equating  the 
derivative  to   0   and  solving  for   K  ,  we  discover 


(7)  K  =  Bn/Dn 


Inserting  this  value  for   K   back  into  (6) ,  we  obtain 


(7)      E  (x,v)  =  A   -  2B2/D   +  B2/D   =  A   -  B2/D 
n         n      n'  n     n   n     n     n   n 


Furthermore,  since  A   does  not  depend  on   x   or   v  , 

2 

minimizing   E  (x,v)   is  the  same  as  minimizing   1  -  B  /(A  D  )  , 
■a    n  -a        n    n  n 


which  is  the  same  as  maximizing 


(8)  R  (x,v)  E  B2/(A  D  ) 

n   '      n   n  n 


Note  that   B    and   D    depend  on   x   and   v   in  (8) ,  although 
the  notation  suppresses  that  fact  to  avoid  being  too  cumber- 
some.  The  factor   A   has  been  retained  because   R  (x,v)   as 
defined  in  (8)  is  just  the  square  of  the  correlation  coeffi- 
cient between  the  two  sequences   C, , . . . ,C    and   d, , . . . ,d   , 

-i         1 '      n        1      n 

which  is  a  valuable  interpretation.   The  LSQ  procedure  selects 

x   and   v   to  maximize   R  (x,v)  ;  i.e.,   x   and   v   are  selected 

n 

to  make  the  theoretical  and  actual  measurements  as  highly 
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correlated  as  possible.   The  quantities   A   ,  B   ,  and  D    are 
r  ^  n     n        n 

all  sums,  so  those  quantities  are  stored  and  updated  at  each 
measurement  (see  Sec.  1.9). 


1.7   Least  Squares  Estimates  (K  unknown) 

Our  strategy  here  will  be  to  eliminate   K   by  sacri- 
ficing the  first  measurement.   To  simplify  the  required  cal- 
culations, assume  that   t,  =  x,  =  0  ;  this  can  always  be  achieved 
by  defining  the  origins  of  time  and  space  to  be  at  the  first 
contact.   From  (1) , 


(9)         (Ci/C1)  1/a   =    |  (x  +  vti  -  xj_)/x|  =  r± 


Ignoring  the  absolute  value  (in  other  words,  guessing  direction) , 
(9)  becomes 


CIO)  x(l  -  ri)  +  vti  -  xi  =  0 


-1/a 
If  the  measured  value   r.  =  (C./C^)       is  substituted  for 

ri   in  (10) ,  the  result  is 


(11)  xCl  -  ri)  +  vtJ.  -  xi  =  ZL    , 


where   €.   is  an  error.   Since  the  left  hand  side  is  a  linear 
i 

...      ?    2 

expression  in   x   and   v  ,  the  task  of  minimizing    i  €. 

i=l 
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with   x   and   v   can  be  accomplished  analytically,  and  we 
define  the  resulting  estimates  to  be  XLS1  and  VLSI, 
respectively. 

A  different  least  squares  estimate  can  be  obtained  if 


(10)  is  divided  through  by   t.   when   i  >_  1  .   The  resulting 

,  and  we  define  the 


sum  of  squares  is    £ 

i=l 


n      1-r .        x . 

+  v  - 


t.         t. 
l         l 


minimizing   x   and   v   to  be   XLS2   and   VLS2,  respectively. 
The  motiviation  for  dividing  through  by   t .   is  that  the  co- 
efficient of  one  of  the  unknowns  (namely  v)  becomes  1.0,  which 
means  that  LS2  estimates  are  even  easier  to  calculate  than  LSI 
estimates.   However,  the  LS2  and  LSI  estimates  are  different, 
and  the  comparisons  in  the  next  section  will  show  that  the  LS2 
procedure  does  not  make  good  use  of  measurements  with  high 
subscripts.   The  reason  for  this  is  that  division  by   t.   gives 
the  high  subscript  terms  very  little  influence  in  the  sum  of 
squares . 

1.8   Comparisons  with   K   unknown 

Figures  17-22  show  a  comparison  of  the  LSI  and  DIR 
procedures  (we  retain  the  name  DIR  even  though   K   is  estimated 
analytically)  in  the  same  circumstances  and  format  as  Figures 
1-16.   Performance  is  of  course  worse  than  in  the  case  where 
K   is  known.   Compare  Figures  2  and  19,  for  example.   For  the 
same  inputs,  the  amount  of  bias  and  especially  the  amount  of 
variance  are  considerably  smaller  in  the  case  where   K   is 
known.   Alternatively,  Figures  2  and  18  exhibit  estimates  of 
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roughly  the  same  quality  where  the  amount  of  noise  is  5  times 
as  large  in  the  case  (Figure  2)  where   K   is  known.   Knowledge 
of   K   is  evidently  worth  a  lot  in  terms  of  the  amount  of  noise 
that  the  system  can  stand. 

The  comparison  between  the  LSI  and  DIR  procedures  leads 
to  the  same  conclusions  as  when   K   is  known.   The  DIR  proce- 
dure seems  to  be  superior  statistically,  particularly  in  the 
low  variability  of  its  estimates. 

Figures  23-28  compare  LS2  with  DIR.   Once  again,  DIR 
wins.   Note  that  the  variance  of  LS2  estimates  declines  slowly 
with  the  number  of  measurements;  the  reason  for  this  was  de- 
scribed earlier  in  1.7.   Compare  Figures  18  and  24.   LSI  and 
LS2  give  the  same  estimates  after  3  measurements,  but  the  LSI 
estimates  improve  much  faster  than  the  LS2  estimates  as  the 
number  of  measurements  increases  beyond  3.   After  20  measure- 
ments, the  LS2 ,  LSI,  and  DIR  procedures  estimate   x   with 
sample  standard  deviations  of  12.64,  3.74,  and  .78,  respectively 
These  are  substantial  differences.   Alternatively,  compare 
Figures  23  and  24  at  the  18 —  measurement.   The  LS2  procedure 
in  Figure  23  has  roughly  the  same  mean  and  variance  as  the  DIR 
procedure  in  Figure  24,  but  the  noise  level  is  10  times  as 
large  in  the  latter.   In  this  instance,  employment  of  DIR  vice 
LS2  is  equivalent  to  a  lOdb  gain  in  sensitivity. 

1.9   An  operational  version  of  Direct  Search 

The  principle  drawback  of  the  Direct  Search  procedure 
is  the  amount  of  time  required  to  make  an  estimate.   In  order 
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to  get  some  idea  of  this  time,  a  version  of  the  Direct  Search 
with  unknown   K   procedure  was  written  in  BASIC  for  the 
TEK  4  052  computer.   The  program  updates  the  quantities   A   , 
B   ,  and   D    used  in  (8)  according  to 


A  ..  =  A   +  E  C2 
n+1    n    n  n 


(12) 


B  _,_,  =  B   +ECd 
n+1    n    n  n  n 

.2 
D  ,n  =  D   +  E  d   ,   where 
n+1    n    n  n 


E   =  exp(-  (t   -  t   , )/t) 
n     ^     n    n-1  ' 


The  reason  for  including  the  factor   E    in  all  of  the  updates 
is  due  to  a  conviction  that  some  sort  of  discounting  mechanism 
should  be  employed  in   practice  to  give  more  weight  to  recent 
than  to  past  measurements,  with   t   being  a  "relaxation  time" 
within  which  the  parameters   x   and   v  might  in  practice  vary 
a  little  bit.   With  this  feature  built  in,  the  program  takes 
about  40  seconds  to  compute,  display,  and  select  the  maximum 
of  the  squared  correlation  coefficient   R    for   50  x  30  =  1500 
grid  points.   Since   R    is  always  between   0   and   1  ,  the  disp] 
is  simply  a   50  x  30  matrix  of  the  first  digit  after  the 
decimal  point  in   R   ;  an  operational  version  might  actually 
draw  a  contour  map.   It  is  felt  that  the  display  provides  valu- 
able information  even  after  the  maximizing  grid  point  is  known, 
since   R   can  have  multiple  local  maxima. 

Figure  29  shows  a  sample  map  after  processessing  3 
observations.   Note  that  the  variables  searched  are  speed  and 
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time  late,  rather  than  speed  and  distance.   The  only  advantage 
of  using  time  instead  of  distance  is  that  the  sign  of  time  is 
necessarily  positive.   The  program  reports  that  the  peak  of 
the  map  is  at  time  late  9.2  and  speed  -11,  which  is  a  point  in 
the  middle  of  the  field  of  9's  in  the  upper  right  hand  corner. 
The  virtue  of  having  the  map  is  that  an  operator  could  see 
that  there  are  also  lots  of  9's  in  the  lower  right  hand  corner; 
i.e.,  it  is  still  possible  that  the  speed  is  positive. 
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TRUE  X=  30.0000  \l=        5.0000  A  L  =   0.0     A2=   0.000L 
XLS1        +/-       VLSI      +/-         XOIR       ♦  /-      VOIR 
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TRUE  X=  30.CC00   V=   5.0000  Al  =   0.0     A2=   0.0010 
XLS1        +/-       VLSI      +/-         XDIR       +/-      VOIR 
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TRUE    X  =     3C.C3C0        \»  =       5.0000    Al=       0.0  A2=       0.0050 

XLS1  *■/-  VLSI  «■/-  XDIR  */-  VOIR  +  /- 
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TRUE    X  =     30.C000       \-       5.0000    Al=       0.010D     A2  =       0.0 
XLS1  +/-  VLSI  +/-  XDIR  +/-  VOIR  +/- 
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TRUE    X=    3Q.C000       V=       5.0000    Al  =       0.0500     \2  =       0.0 
XLS1  +/-  VLSI  +/-  XDIR  +/-  VOIR 
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TRUE    X=    30.CC00       V=       5.0000    Al=       0.0100    A2=       0.0010 
XLS1  +/-  VLSI  +/-  XDIR  ♦  /-  VOIR 
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TKLE  X=     30.  C  CCO       \i=       5.00  0  0    Al  =       0.0             A2=       0.30  01 

XLS2  +  /-                VLS2  +/-                     XDIR  ♦•/-  VOIR                   +/- 
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TRUE    X=     3C.CC00        v=        5.0000    Al=       0.0  A2=       0.0010 

XLS2  +/-  VLS2  +/-  XDIR  */-  VOIR 
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TRUE    X=     30.0C00        V=       5.0000    41=      0.0  A2  =       3.0050 

XLS2  ♦/-  VLS2  +/-  XDIR  +/-  VOIR 
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2.    Maximum  likelihood  estimate  for  multiplicative  noise. 

2.1  Introduction 

As  was  pointed  out  in  Section  1.2  the  estimate  obtained 
by  the  direct  search  method  has  the  virtue  of  being  a  maximum 
likelihood  estimate  in  case  of  an  additive  white  Gaussian  noise. 
If  a  multiplicative  noise  is  present,  however,  this  is  no  longer 
the  case.   It  is  therefore  of  some  interest  to  investigate  to 
what  extent  the  estimation  may  be  improved  by  using  the  true 
maximum  likelihood  estimator.   To  do  so  it  is  of  course  neces- 
sary to  make  an  assumption  about  the  distribution  of  the  multi- 
plicative noise.   We  shall  assume  that  the  multiplicative  noise 
is  log-normal  just  as  we  have  done  in  Section  1.   Although  at 
present  there  seems  to  be  only  a  meager  evidence  in  favor  of 
this  assumption,  there  is  not  much  evidence  in  the  contrary 
either.   Thus,  at  present,  our  assumption  can  be  perhaps  justified 
as  a  convenient  working  hypothesis. 

2.2  The  model 

We  assume,  as  in  Section  1,  that  the  true  range  to  the 
target  at  the  time  of  the  i-th  measurement  is 


CD  ri  =  |u  +  vti  -  xi|  ,  i  =  l,...,n  , 


where   v   and   u   are  the  speed  and  initial  distance  to  a  target. 

The  position   x.   of  the  sensor  at  the  time   t.   of  the  i-th 
r          l                              i 

measurement  is  assumed  to  be  measured  with  respect  to  a  reference 
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direction  along  the  target  track  with,  the  initial  position  at 
t,  =  0   of  the  sensor   x,  =  0   as  the  origin. 

The  target  is  assumed  to  travel  with  a  constant  speed 
and  not  be  overtaken  by  the  sensor  during  the   n   measurements. 
Thus,  if  the  target  travels  in  the  reference  direction  the 
quantities   u,v,  t.   and  u  +  vt.  -  x.   in  (1)  are  all  positive, 
while  if  the  target  travels  against  the  reference  direction  the 
same  effect  is  accomplished  by  changing  the  signs  of  all   x.  *s  . 
Hence  if  we  define  for  all   i  =  l,...,n 


xi  =  i 


x.     if  the  target  travels  in  the  reference 
direction, 


-x .     if  the  target  travels  agains  the 
reference  direction, 


we  can  write  for  the  range 


(2)  r .  =  u  +  vt.  -  x.  ,  i  =  l,...,n  , 


with   r.  ,  u  ,  v  ,  t.   all  positive  in  either  case.   This  elimi- 
nates the  cumbersome  absolute  value  in  (1) . 

The  i-th  sensor  measurement   C.   is  assumed  related  to 


the  range   r.   by  the  equation 


(3) 


Ci  =  Kri~aQi  '  i  =  1' *  *  * 'n  ' 
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where  a    >  0   is  known  constant,   K  >  0   is  unknown  constant 
signal  strength,  and   Q-,,...,Q    are  independent,  identically 
distributed  positive  random   variables  representing  multipli- 
cative measurement  noise.   The  distribution  of   Q .   is  assumed 


log-normal,  i.e 


aW 


Q.  =  e 

1 


where   a   is  a  constant  and  W.   are  standard  Gaussian  random 
variables.   No  additive  noise  is  considered  in  this  model. 


2.3   Derivation  of  the  maximum  likelihood  estimator 

Let  the  random  variables   Q.   in  (3)  have  the  log-normal 


density 


(4) 


fQ(q)  =  ~7F  5  e 


1   o  2 
— 7?   In    q 


,  q  >  0  , 


where   a  >  0   is  a  parameter. 

AlWi         2 
Remark ;   This  is  equivalent  to   Q.  =  e      with   a 


2    2 

A    A 

e  1(e  x  -  1) 


and  W.   standard  normal, 
l 

The  likelihood  function  is  from  (3)  and  (4) 


(5) 


I  V«1  /  •  •  •  /^-'T_./ 


n  r . 

=  rr 


i=l 


K 


C.r. 

i  i 


K 


n 


TT    1    1 

IT  ——  q-   exp 
i=l  a/2TT   i 


2a' 


In 


-    C  .  r  .  a 
2      i  l 

K 
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and  this  is  to  be  maximized  with,  respect  to   u  ,  v  or  u  ,  v  ,  K 
for   K   unknown  subject  to  the  constraints 


r .  =  u  +  vt . 

1  1 


x.  ,  i 


1 ,  .  .  .  ,  n 


Taking    logarithm  we   get 


n 


(6) 


Znf(C,,...,C    )    =  -  n   In   a/Tn   -      J      Jin  C . 
1  i=l 


-   -i      I      In 
2a      i=l 


n  0    C. r . 

2      li 


a 


K 


so    that    the   problem   is    equivalent   to  minimizing 


(7) 


n  -    C.r . 

v       o~2      i    i 


i  =  l 


K 


If  the  signal  strength   K   is  known  we  thus  conclude 
that  the  maximum  likelihood  estimate   u  ,  v  minimizes  the 
objective  function 


n 


(8) 


D(u,v)  =  I 
i=l 


in 


C.r  . 
l  l 


a 


K 


subject    to    the    constraints 


r.    =   u   +   vt.    -   x.     ,    i    =    l.....n 
l  li  '         ' 
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For  the  case  of  unknown   K   we  note  that  the  function 


to  be  minimized 


n     ~  c  .  r  .      n 
I      In2   -^  =  I 


i=l 


K 


i=l 


in   c.r.   -  £nK 

1  1 


has  as  a  function  of   K  >  0   a  single  minimum  at 


K  = 


[    n                  1 

1 

TT  Cr.a 
.!  ',      11 

n 

1=1 

Upon   substitution   back    into    (7)    we   obtain    the    objective    function 


n 


(9) 


D(u,v)     =      I 
i=l 


Jin   C.r.a 
l    l 


where  we  have  denoted 


r  .  =  r  . 
i     i 


n 


*""  n 


C. 

l 


c . 

1 


n 


TT  c 


i 

n 


i  =  1, . .  .  ,  n 


Thus  the  maximum  likelihood  estimate   u,  v   for  the  case  of 
unknown  signal  strength  K   must  minimize  (9)  again  subject  to 
the  constraints 


r  .  =  u  +  vt . 
i  l 


x.  ,  i 


=  1 ,  .  .  .  ,  n 
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Notice  that  in  either  case  the  estimator  is  still  of  a 
least-squares  type  but  the  dependence  on  the  parameters   u,  v 
is  nonlinear.   Hence,  for  each  set  of  observations   C,,...,C 
the  maximum  likelihood  estimate   u,  v  must  be  determined  by 
numerical  minimization.   In  view  of  the  nature  of  the  problem 
considered  in  this  report  the  simplest  approach  is  to  use  the 
direct  search  method  as  in  Section  1,  applied  this  time  to  the 
objective  function  (8)  or  (9). 

2.4   Results  of  computer  simulation 

Since  the  purpose  was  to  investigate  the  possible  im- 
provement due  to  the  use  of  maximum  likelihood  estimator  simu- 
lation was  performed  with  the  same  set  of  parameters  as  in 
Section  1.4  (for  known   K)  and  Section  1.8  (for  unknown   K) . 
The  data  are  displayed  in  the  same  format  in  Figures  30  to  35. 
For  the  case  of  known   K  ,  Figure  30  should  be  compared  with 
Figure  9,  Figure  31  with  Figure  11,  Figure  32  with  Figure  12. 

It  is  seen  that  indeed  the  maximum  likelihood  estimates 
result  is  smaller  standard  deviations  and,  especially  for  larger 
noise  intensity,  less  bias. 

Figure  33  was  included  to  see  how  much  would  the  per- 
formance deteriorate  if  the  observation  noise  were  in  fact  ad- 
ditive rather  than  multiplicative.   Comparison  with  the  XDIR 
and  VDIR  columns  of  Figure  5  shows  that  one  must  pay  an  esti- 
mation penalty  for  assuming  multiplicative  noise  when  the  noise 
is  actually  additive.   Nevertheless,  the  maximum  likelihood 
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estimator  for  the  multiplicative  noise  still  performs  better 
than  the  least  squares  algorithm. 

For  the  case  of  unknown  signal  strength   K   the  results 
of  simulation  are  in  Figures  34  and  35  which  should  be  compared 
with  Figures  20  and  21  respectively.   Here  the  performance  of 
the  maximum  likelihood  estimator  is  about  the  same  as  the  direct 
search  method  although  somewhat  smaller  bias  can  be  detected  in 
Figure  35. 

In  conclusion,  if  the  noise  were  indeed  multiplicative 
log-normal  and  if  no   additive  noise  were  present  then  the  maxi- 
mum likelihood  estimator  would  be  generally  preferable.   However, 
the  direct  search  method  seems  to  be  more  robust  and  would  prob- 
ably be  safer  to  use  at  least  until  a  verifiable  distributional 
assumptions  can  be  made  about  the  nature  of  the  noises  involved. 
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TRUE    X=     30.0000       V=       5.0000    Al  =       0.1000    A2  =       0.0 


XDIR 


+  /- 


VOIR 


+  /- 


30.  26 

2.09 

3.45 

1  .57 

29.  64 

2.00 

7.18 

6.61 

29.  54 

2.0  3 
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5.95 
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2.0  9 

4.99 

4.57 

29.  98 

2.20 

5.11 

3.49 

30.  13 

2.37 

4.70 

3.28 

30.  06 

2.19 

4.98 

2.70 

30.08 

1.90 

4.  9  9 

2.11 

30.06 

1.94 

4.99 

1.32 

30.22 

1.62 

4.77 

1.36 

29.96 

1.61 

5.10 

1.26 

29.  96 

1.43 

5.12 

0.96 

30.02 

1.32 

5  .  02 

0.83 

30.  10 

1.32 

4.97 

0.81 

30.  04 

1.26 

4.96 

0.63 

30.06 

1.21 

5.02 

0.56 

30.  18 

1.29 

4.95 

0.57 

30.  12 

1.28 

4.95 

0.52 

30.22 

1.42 

4.91 

0.55 

30.  18 

1.16 

4.94 

0.43 

Figure    30 


TRUE    X=     30.C000       V=        5.0000    Al=       0.5000    A2=       0.0 


XDIR 


♦  /- 


VOIR 


+  /- 


31.44 

9.78 

4.0  2 

3.23 

31.  18 

3.73 

7.96 

7.05 

31.  12 

8.09 

8.07 

7.05 

30.  20 

7.54 
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6.66 
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7.09 

5.79 
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30.08 
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30.08 
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30.  70 
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2.63 

30.  50 

5.75 

A. 33 

2.50 

30.68 

5.61 
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2.19 

30.48 

5.35 

4.90 

1.96 

Figure    31 
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TRUE  X=  30.0000   V=   5.0000  Al=   1.0000  A2=   0.0 


XDIR       +  /-      VOIR       ♦  /- 


31.  80 

14.43 

5.33 

5.06 

32.36 

13. 5o 

3.46 

7.14 
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3.43 

Figure  32 


TRUE  X=  30.0000   V=   5.0000  Al=   0.0     A2=   0.0200 


XDIR       +/-      VOIR       +/- 


32.  54 

12.23 
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Figure  33 
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TRUE    X=    30.0000       V=       5.0000    AL=      0.0100    A2  =       0.0 


XOIR 


+  /- 


VOIR 


+/- 


50.00 

0.04             13.72 

2.36 
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Figure    34 

0.10 

TRUE    X=    30.0000       \,=       5.0000    Al  =      0.0500 


A2=       O.Q 


XDIR 


+  /- 


VOIR 


♦  /- 
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Figure   35 
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3.    Determination  of  target  direction. 

3.1  Introduction 

For  the  purpose  of  target  tracking  it  is  important  to 
determine  the  direction  of  target  motion  after  only  few  initial 
measurements.   This  should  be  done  with  as  little  knowledge  of 
the  signal  parameters  as  possible.   Specifically,  the  target 
speed,  distance  and  the  signal  strength  will  typically  be  un- 
known at  the  beginning  and  cannot  be  reliably  estimated  from 
only  a  few  measurements.   In  this  section  we  propose  a  simple 
method  for  deciding  the  target  direction.   The  results,  however, 
are  only  preliminary.   More  investigation  is  needed  to  develop 
the  method  into  a  workable  algorithm. 

3.2  Description  of  the  method 

We  assume  here  that  the  direction  of  the  target  is  to 
be  determined  immediately  after  three  initial  measurements.   We 
also  assume  that  the  locations  of  the  two  measurements  follow- 
ing the  initial  one  bracket  the  latter  and  that  the  coordinate 
axis  is  oriented  towards  the  last  measurement.  (See  Figure  36.) 


FIRST 
MEASUREMENT 


INITIAL 
MEASUREMENT 


SECOND 

MEASUREMENT 


COORDINATE  AXIS 


-0 


(REFERENCE   DIRECTION) 


Figure    36 
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Thus  the  initial  measurement  is  made  at  the  coordinate 

xn  =  0  ,  the  first  at   -  x,  <  0   and  the  second  at   x?  >  0  . 

Note  that  x,   and   x~   are  distances  rather  than  coordinates, 

hence  both  are  positive.   The  times  are   t-  =  0  ,  t,  >  0   and 

t~  >  t,   respectively.   Let  the  initial  distance  to  the  target 

be   u  >  0   and  the  target  speed   v  >  0  ,  let   r.  ,  j  =0,  1,  2, 

be  the  distance  to  the  target  at  the  time   t .  .   We  then  have 

1 


rQ  =  u 

r,  =  u  +  vt,  +  x, 

r2  =  u  4-  vt2  -  x2 


if  the  target  travels  in  the  reference  direction,  and 


r0  =  U 

r,  =  u  +  vt,  -  x, 

r2  =  u  +  vt2  +  x2 


if  the  target  travels  against  the  reference  direction, 
If  we  look  at  the  second  difference 


A2rQ  =  Cr2  -  rx)  -  Lr±   -  rQ) 


we  see  that  it  can  be  written  as 


A2rQ  =  v(t2  -  2t1)  +  Cx2  +  2xx) 
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Thus  if   t2  =  2t,   the  second  difference  is  positive  or 
negative  depending  on  whether  the  target  travels  against  or 
along  the  reference  direction  respectively.   Note  that  it  is  in 
fact  sufficient  that  the  term   v(t~  -  2t,)   be  small  in  magni- 
tude relative  to   x~  +  2x,  ,  a  condition  which  can  always  be 
accomplished  since  the  times   t,  ,  t~   and  distances   x,  ,  x? 
depend  only  on  the  motion  of  the  observation  platform. 

Now  the  measurements   C~  ,  C, ,  C_   are  related  to  the 
distances  to  the  target  by 


C .  =  Kr  .  aQ  .  +  W  . 

3    333 


J  =  0,1,2, 


where   K   is  the  signal  strength  and   Q.'s   and  W.'s   are  random 

variables  representing  the  measurement  noise.   If  the  noise  were 

absent  CQ  ■  =  1   and  W.  =0)  we  could  solve  for   r.   and  calculate 
3  3  3 

the  second  difference 


A  rQ  =  K 


-  2C 


+  C, 


1 


The  direction  of  the  target  motion  would  then  be  deter- 
mined from  the  sign  of  the  expression 


T  =  C, 


-  2C 


1 

a 


+  a 


In  particular  if   T  <  0   the  target  travels  in  the 
reference  direction  while  if   T  >  0   the  target  travels  against 
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the  reference  direction.   Notice  that  neither  the  signal  strength 
K   nor  the  true  target  speed   v   needs  to  be  known  to  determine 

the  direction,  only  the  value  of  the  exponent  a      is  required. 

2 
In  the  presence  of  noise  the  second  difference   A  rQ 

is  no  longer  proportional  to  the  test  statistic   T   but  we  can 

still  hope  that   T   will  do  the  job.   In  fact,  if   W.  =  0   and 


LQ3 


1/a 


=  1   then   T   is  a  point  estimator  for   K     A  r. 


obtained  by  the  method  of  moments.   If  we  make  a  natural  assump- 
tion that  the  mean  of  the  multiplicative  noise  equals   1   and 


if  the  noise   Q  .   is  log-normal  then   E 

] 


A 


1/a 


=  (E 


LQDJ 


l/a: 


is  indeed  one. 

The  next  section  describes  the  results  of  a  simulation 
experiment  designed  to  test  the  usefulness  of  the  statistic   T  . 

3.3   Computer  simulation 

In  order  to  gain  some  appreciation  of  the  performance 
of  the  method  we  again  used  computer  simulation.   The  parameters 
were  the  same  as  in  Section  1,  true  range  30NM,  target  speed  5KT, 
signal  strength   K  =  1   and  the  exponent  a   =  1  .   For  each 
simulation  run  1000  samples  of  the  first  three  observations 
C0,  C,,  C2   were  generated.   Each  sample  was  generated  both  for 
positive  true  target  direction  (along  the  reference  axis)  and 
for  negative  true  direction  (against  the  reference  axis.)   Errors 
in  direction  determination  were  counted  and  the  counts  were  used 
to  estimate  the  error  probabilities  for  either  of  the  true 
directions . 
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For  the  first  five  runs  we  used  the  noise  model  of 
Section  1,  equation  (5),  i.e.  multiplicative  log-normal  noise 
and  additive  Gaussian  noise.   This  was  done  mainly  for  the  sake 
of  comparison  with  earlier  simulation  results  in  this  report. 
The  results  are  summarized  in  TABLE  1.   It  is  seen  that  for 
moderate  values  of  the  coefficients   A,   and   A-   the  error 
probabilities  are  small  enough  for  the  method  to  be  useful.   For 
larger   A,   and   A2  ,  however,  the  direction  determination  tends 
to  be  random,  which  is  not  too  surprising  with  only  three 
observations . 

For  the  next  series  of  simulation  runs  the  observations 
were  generated  using  the  equation 


C..    =  KC1  +  Zj>rj~a  »     j  =  0,  1,  2, 


where   Z_,  Z,,  Z ~   were  independent  random  variables  with  the 
density 


fz(z)  =  ±^-l    (1  -  |2|B)  ,     |z|  <  1 


The  parameter   3  >  0   was  varied  from  about  300  to  zero.   Note 
that  the  density  function  is  symmetric  about  zero  and  as   3 
decreases  it  flattens  towards  the  uniform  density  on  (-1,  1) . 
Thus  the  multiplicative  noise  factor   1  +  Z .   always  averages 
to   1   but,  as   3   decreases,  tends  to  distort  each  true  meas- 
urement Kr.~a   by  larger  and  larger  percentage,  rendering  it 
almost  useless  for   3  =  0  (.uniform  distribution)  . 
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The  estimated  error  probabilities  are  shown  in  Figure  37, 
On  the  abscissa  of  the  graph  is  the  ratio   of  the  noise  standard 
deviation   a(l  +  Z.)   to  that  of  the  uniform  distribution  on 
(-1,  1) ,  viz.   /3  ,  which  seems  to  be  more  informative  than  the 
parameter   g   itself.   Of  course,   g   can  be  easily  recovered 
since 


°  U  +  V  -  (6  +  2)  (g  +  3) 


TRUE  DIRECTION: 

POSITIVE 

NEGATIVE 

Al 

A2 

ERROR  PROBABILITY 

.001 

.001 

.00 

.00 

.01 

.001 

.01 

.00 

.05 

.001 

.05 

.06 

.1 

.001 

.27 

.22 

.5 

.001 

.47 

.43 

TABLE  1. 
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Notice  that  for  the  problem  parameters  used  here  there 
is  a  bias  towards  the  positive  direction.   However,  the  results 
do  indicate   that  the  method  may  lead  to  a  practical  direction 
determination  algorithm. 
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4.    Equivalent  sweep  width. 

In  this  section  we  consider  the  two  dimensional  problem 
of  detecting  a  target  that  will  be  detected  if  either 

a)  the  searcher  comes  within   r   of  the  target  in  any 
direction. 

b)  the  searcher  crosses  a  line  segment  of  length   L 
that  the  target  trails  behind  it,  where   L  >  r  . 

The  situation  is  illustrated  in  Figure  1.   The  angle   9   is  the 
angle  perpendicular  to  the  searcher's  course  in  target  centered 
coordinates.   The  heavy  line  of  length  w(8)   is  the  sweep  width 
in  the  sense  that  detection  will  occur  if  and  only  if  the 
searcher's  track  crosses  it.   The  formula  for   w(6)   is 


(1) 


w(6)  =  { 


r   +   Lsint 


2r 


•  \ 


>  90  e  arcsml- 


if  e  < 


There  are  two  assumptions  commonly  made  about   8  :   one 
is  that   9   is  uniformly  random  (area  search) ,  and  the  other 
is  that   9   is  determined  by  relative  speed  in  a  barrier  search. 
These  two  cases  are  analyzed  below. 

4.1   Area  search 

If   9   is  completely  random,  the  average  value  of   w(9) 
can  be  found  by  assuming   9   to  be  uniformly  random  in  the  first 
quadrant : 
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w  = 


U0  7T/2 

/  2rd9  +   /   (r  +  Lsin6) d( 

o       eQ 


2r9Q  +  r(iT/2 


i~)  +  Lcos90 


(2) 


7T 


+  AF^7- 


when   r  =  0  ,  the  equivalent  sweep  width  is  (2/tt)L  ,  a  familiar 
result.   The  point  is  that  there  is  also  an  equivalent  sweep 
width  in  the  case  where   r  ^  0  . 


4.2   Barrier  Search  (back  and  forth  perpendicular  to  target  tracks 

/~2  2 

In  this  case   8   is  arcsin  v//v  +  u   ,  where   u   and  v 

are  the  target  and  searcher  speeds,  respectively.   Let  the 
result  of  substituting  this  into  (1)  be   W   ,  since  we  wish  to 
emphasize  the  effect  of  searcher  speed: 

f  n       2 

+  Lv//v   +  u      i 


(3)     W   =  { 

v 


r  +  Lv/i 


2r 


f    r/L  <  v/vV2  +  u2 

n       2 

if     r/L  >  v//v   +  u 


The  probability  of  detection  when  the  searcher  must  cover  a 
length  of  barrier   B   is  [reference  1] : 


w 


(4) 


P  .,  =  mi 


i    v  /  ~  2  .  2 
n  \1,    -5-  /l  +  v  /u  - 


3 
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For  example,  Figure  2  shows  a  graph  of   P,   vs.   v/u   in  the 
case   r=l,L=2,B=10/  obtained  by  substituting  (3)  into 
(4)  .   The  point  is  that   P  ,   is  a  very  strong  function  of   v/u 
over  the  region  where   P-,  <  1  —  the  curve  is  actually  convex. 
This  strong  dependence  on  speed  is  because  sweep  width  increases 
with  speed  due  to  an  increasingly  favorable  geometry,  so  the 
all  important  product  of   v   and  W    increases  fast.   In  prac- 
tice, increasing  speed  will  also  have  a  detrimental  effect  on 
signal-to-noise  ratio,  so  the  "always  proceed  at  top  speed" 
conclusion  indicated  by  this  simple  analysis  does  not  hold. 
Nonetheless,  the  geometrical  considerations  are  important.   We 
claim,  for  example,  that  the  best  speed  in  a  barrier  patrol  is 
greater  than  the  best  speed  in  area  search,  since   W   in  (2) 
does  not  depend  on   V. 
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